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, Abstract. The POINT- AGAPE collaboration is carrying out a search for gravitational microlensing toward M31 

to reveal galactic dark matter in the form of MACHOs (Massive Astrophysical Compact Halo Objects) in the 
halos of the Milky Way and M31. A high-threshold analysis of 3 years of data yields 6 bright, short-duration 
microlensing events, which are confronted to a simulation of the observations and the analysis. The observed 
signal is much larger than expected from self lensing alone and we conclude, at the 95% confidence level, that at 
least 20% of the halo mass in the direction of M31 must be in the form of MACHOs if their average mass lies in 
the range 0.5-1 M©. This lower bound drops to 8% for MACHOs with masses ~ 0.01 Mq. In addition, we discuss 
a likely binary microlensing candidate with caustic crossing. Its location, some 32' away from the centre of M31, 
supports our conclusion that we are detecting a MACHO signal in the direction of M31. 
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1. Introduction by a factor ~ 5 to be accounted for by known popula- 

■ 1 tions of stars . Indeed, further analysis recently confirmed 

Gravitational microlensing, as first noted by | Paczynski| these results feennett et all200dlBejmetll20o3l . This ex- 

(1986D, is a powerful tool for the detection of massive cesg ig attributed to MACHOs of mass ~ 0.4 M in the 

astrophysical halo compact objects (MACHOs), a pos- MW halo contributing f „ 2 0%, altho ugh this resu lt has 

sible component of dark matter halos. Observations to- been challenged by se veral authors fe.g lJetzer et ahlEool 

ward the Magellanic Clouds by the first generation of mi- | B elokurov et allEioi . These exciting and somewhat con- 

crolensing surveys yielded important constraints on the tradictory results challenge us to probe the MACHO dis- 

Milky Way (MW) halo. The EROS collaboration ob- tribution along different MW lines of sight and in different 

tained an upper limit (/ < 20%) on the contribution b y Ea l ax ies 
MACHOs to a standard MW halo l)Afonso et alJl2003jh 

and the results of their latest analysis strengthen this con- M31, being both nearby and similar to the MW 



elusion l|Tisserand fc Milsztainll2005l). Also, ac cording to is a suitable targ et for such a search 1 Crottsl Il992t 

the MACHO collaboration l|Alcock et al.ll2000j) . the opti- iBaillon et al"]ll993|h It allows us to explore the MW halo 

cal depth toward the Large Magellanic Cloud is too large along a different line of sight. It has its own halo that 

can be studied globally, and its high inclination is ex- 

* UMR 7164(CNRS, Universite Paris 7, CEA, Observatoire pected to give a strong gr adient in the spatial dist ribu- 

de Paris) tion of microlensing events l|Crottslll992t IJetzerl l994'l . We 
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note, however, that the latter feature, which was at first 
believed to provide an unmistakable signature for M31 
microlensing halo events, seems to be shared, at least to 
some extent, by the variable star population within M31 
l|An et alJl2004a|) . 

Several collaborations have und ertaken searches fo r 
microlensing towar d M31: AGAPE llAnsari et alJ ll9991. 
SLOTT-AGAPE lEaTchi Novati et all l2003h. MEGA 
dde Jong et alJ l2004h. Columbia- VATT ()Uglesich et alJ 
\2004. WeCAPP llRiffeser et alJ l2003h and NainiTal 
ll.Toshi et all 1200,1) . Up to now, while some microlens- 
ing events have been detected, no firm conclusion 
about their physical meaning has been reported. In 
particular, the POINT-AGAPE collaboration has pre- 
sented a first analysis focused on t he search for bright, 
short-duration microlensing: e vents l)Auriere et alJ EfJOll 
IPaulin-Henriksson et al.ll2003h . 

In this paper, we report the first constraints on the 
MACHO fraction of the combined MW and M31 halos 
along the line of sight to M31. We give a complete account 
of our systematic search for bright short-duration events, 
present the 6 selected microlensing events, and then de- 
scribe the simulation used to predict the characteristics 
of the expected events and their frequency. We proceed 
in two steps: a Monte Carlo simulation produces an ini- 
tial (quantifiably over-optimistic) estimate of the number 
of expected events, then a simulation of events (hereafter 
referred to as "event simulation") on the actual images 
allows us to assess the detection efficiency of the analy- 
sis pipeline for the type of events produced by the Monte 
Carlo. 

In the search for a MACHO signal we must deal with 
two main backgrounds: (i) variable stars masquerading as 
microlensing events and (ii) self-lensing events (for which 
both the lens and the source are part of the luminous 
components of M31 or MW). We eliminate the first (see 
below) and partially isolate the second using their distinc- 
tive spatial distribution. 

The paper is organised as follows. In Sect. |3 we re- 
call the observational setup and then describe our analysis 
pipeline. The detected microlensing signal is discussed in 
Sect. 13 In Sect. 01 we describe the Monte Carlo simulation 
of the experiment and describe its predictions. In Sect. 
03 we evaluate the detection efficiency of the pipeline. In 
Sect. we summarise the analysis and discuss what con- 
clusions can be drawn about the fraction of M31 and MW 
halos in the form of MACHOs. 



2. Data analysis 

2.1. Setup, data acquisition and reduction 

In this work we analyse data acquired during three sea- 
sons of observation using the Wide Field Camera (WFC) 
mounted on the 2.5m Isaac Newto n Telescope (INT) 
l|Auriere et alJ l200li lAn et all 120044 - A fourth year of 
data is currently being analysed. Two fields, each ~ 0.3 
deg 2 , north and south of the M31 centre are monitored 




Fig. 1. Projected on M31, we display the boundaries 
of the observed fields (red lines), and the centre of M31 
(cross). Circles mark the positions of the 6 microlensing 
events issued from the selection pipeline (Sect. I3~T1) . The 
open circle (S4) corresponds to the event seen toward M32. 
The star (S5) to the binary event candidate discussed in 
Sect. ESI 

(Fig. The data are taken in two passbands (Sloan r 
and either Sloan g or Sloan i), with exposure time be- 
tween 5 and 10 minutes per night, field and filter. Each 
season of observation lasts about six months, but with 
very irregular sampling (especially during the third one). 
Overall, for r data, we have about 120 nights of observa- 
tion. At least two exposures per field and filter were made 
each night with a slight dithering. Although they are com- 
bined in the light curve analysis, they allow us to assess, 
if necessary, the reality of detected variations by direct 
inspection of single images. 

Data re ductio n is performed fol lowing 
Ansari et all lll997h. fCalchi Novati et all l|2002h and 
Paulin-Henriksson et alJ <|2003j) . Each image is geometri- 
cally and photometrically aligned relative to a reference 
image (one per CCD, the geometric reference being the 
same for all the filters). Ultimately, in order to deal with 
seeing variations, we substitute for the flux of each pixel 
that of the corresponding 7-pixel square "superpixel" 
centred on it, the pixel size being 0.33", and we then 
apply an empirical correction, again calibrating each 
image against the reference image. 



2.2. Analysis: selection of microlensing events 

To search for microl ensing events , we u s e the " pixel- 
lensing" techniqu e (|Baillon et all Il993t iGouldl Il99fit 
lAnsari et a,l.lll99^ . in which one monitors the flux vari- 
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ations of unresolved sources of each pixel element of the 
image. 

A difficulty, specific to pixel lensing, is that genuine 
microlensing events might be polluted by one of the nu- 
merous variable objects present in the neighbouring pixels. 
To avoid loosing too many bona fide microlensing events 
while accounting for the variable background, we look for 
microlensing-like variations even on those light curves on 
which a second bump is detected. In particular, on each 
light curve we first look for and characterise mono-bump 
variations for each season separately. Only as a final step 
do we test for bump uniqueness on the complete light 
curve in a loose way as explained below. This test allows 
for the presence of variable stars within the superpixel con- 
taining the lensed source and so, as a bonus, in principle 
could allow us to detect microlensing of variable objects. 

In addition to the physical background of variable 
stars, the search for microlensing-like flux variations, in 
particular the short ones, is plagued by the detection of 
"fake" variations, mainly due to bad images, defects on the 
CCD, saturated pixels associated with extremely bright 
stars, a nd cosmic rays (thes e issues are discussed in more 
detail in lTsanras et alfeOOfll) . The only safe way to remove 
these artefacts is to visually inspect the images around 
the time of maximum, although there may be other use- 
ful hints, such as an anomalous distribution of the times 
of maximum or in the spatial distribution. To obtain a 
"clean" set of variations we first run the complete pipeline, 
identify and remove bad images, and mask bad pixels. 
Then, we rerun the analysis from scratch. 

Before proceeding further with the pipeline, we mask 
a small region right around the centre of M31, ~ 1' x 
1', where, in addition to problems caused by saturation, 
the severe uncertainty in modelling the experiment would 
prevent us from drawing any significant conclusion about 
the physical implication of any result we might obtain. 

As a first step, we establish a catalogue of signifi- 
cant flux variations (using the r band data only, which 
are both better sampled and less seriously contaminated 
by intrin sically variable star s than the i band data). 
Following ICalchi Novati et alJ l)2003h . we use the two es- 
timators, L and Q, which are both monotonic functions 
of the significance of a flux variation, to select candi- 
dates. Note t hat the previous POINT- A G APE selections 
present e d in IPaulin-Henriksson et alJ l)2003(k lAn et all 
I 2004a JkrBelokurov et a.1.1 (|2005^ . have been carried out us- 
ing the L estimator only. 

We define 

L = - In (IL jelmmp P(®\$ >$j)) given $ hfcs , a 3 , (1) 



where 




&j and <jj are the flux and associated error in a super- 
pixel at time tj, &bkg is an estimator of the baseline level, 
defined as the minimum value of a sliding average over 18 



epochs. A "bump" is defined as a positive variation with at 
least 3 consecutive points more than 3cr above the base- 
line, and it is regarded as ending after two consecutive 
points less than this threshold. We define 

2 9 
V — V 

q svconst ^pacz /,^\ 

where xlonst ^ s calculated with respect to the constant-flux 
hypothesis and x 2 pacz i s the x 2 calculated with respect to a 
Paczynski fit. Let us stress that Q is evaluated for each full 
season, while L is evaluated only inside the bump. At this 
point, we keep only light curves with Q > 100. Since Q 
is biased toward mono-bump variations, this step allows 
us to remove the unwanted background of short-period 
variable stars. 

Although it has already been described in 

ICalchi Novati et all 120021 , we return to a crucial 
step of the above analysis. For each physical variation, 
there appears a whole cluster of pixels with Q > 100 
(with typical size range from 4 to 30 pixels). From the 
Q values of all light curves, we construct a Q map for 
each season. We then proceed to the actual localisation 
of the physical variations 1 . First we identify the clusters 
(which appear as hills on the map). Then we locate 
the centre of the cluster as the pixel with the highest 
value of the L estimator. The main difficulty arises 
from the overlap of clusters. Indeed we must balance 
the search for faint variations with the need to separate 
neighbouring clusters. In the following, we will refer to 
this crucial part of the analysis as "cluster detection". 
It must be emphasised that this step cannot be carried 
out on separate light curves, but requires using Q maps. 
The impossibility of including this cluster detection in 
the Monte Carlo (Sect. gives us one of the strongest 
motivations for the detection efficiency analysis described 
in Sect. |SJ After the clusterisation, we are left with ~ 10 5 
variations. 

The following part of the analysis is carried out work- 
ing only on pixel light curves. 

As a second cut, we remove flux variations having too 
small a signal-to-noise ratio (most likely due to noise) by 
demanding L\ > 40, L\ being associated with the bump. 
If the light curve shows a second bump over the three 
seasons, characterised by L2, we then demand that this 
satisfies L2 < 0.5Li. As we are only looking for bright 
bumps (see below), we consider such a significant second 
bump to indicate that these bumps most likely belong to 
a variable star. 

We estimate the probability for the lightcurve of a 
given event to be contaminated by a nearby variable 
source as the fraction of pixels showing a significant vari- 
ation, L\ > 40. This fraction stronlgy depends on the dis- 
tance from the centre of M31: from ~ 10% — 20% in the 
inner M31 region, within an angular radius of 8', down to 
~ 8% in the outer region. 

1 We use here a software developed within the AGAPE col- 
laboration. 
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We characterise the shape of th e varia tion by studying 
its compatibility with a IPaczvhskll l|l98rj) shape. We per- 
form a two-band 7-parameter fit: the Einstein time, tg, the 
impact parameter, Uq, the time at maximum magnification 
to, and the band dependent flux of the unresolved source, 
<fi* , and the background flux, (fib, of the bump in each of 
2 bands (r and either i or g according to the available 
data along the bump) 2 . Throughout the analysis we use, 
as an observable time width, the full-width-half-maximum 
(FWHM) in time of the Paczyhski curve, fi/2, and the flux 
increase A<£> of the bump, both of which ar e functions o f 
the degenerate parameters tg, uq and <fi* (|Gouldl fl996) . 
Using the flux deviation in the two bands, we evaluate in 
the standard Johnson/Cousins magnitude system i?(A$), 
the "magnitude at maximum" of the bump, and its colour, 
either V — R or R — I. The simultaneous Paczyhski fit in 
two bands effectively provides a test of the achromaticity 
expected for microlensing events. 

As a third cut, we use the goodness of the Paczyhski 
fit as measured by the reduced x 2 - For short events, the 
behaviour of the baseline would dominate the \ 2 ■ To avoid 
this bias, we perform the fit in a smaller "bump region" de- 
fined as follows. A first Paczyhski fit on the whole baseline 
provides us with the value of the baseline flux (fib and first 
estimates of the time of maximum magnification to and 
the time width ti/2- Using these values we compare two 
possible definitions of the bump region and use whichever 
is the larger of: (i) the time interval inside to ±3ti/ 2 , and 
(ii) the time interval that begins and ends with the first 
two consecutive points less than 3 a above the background 
on both sides of to- The final Paczyhski fit is carried out 
in this "bump region" with the basis flux (fib fixed in both 
colours, and this fit provides the values of the 5 remaining 
parameters. 

Our third selection criterion excludes light curves with 
X 2 /dof > 10. 

We fix this threshold high enough to accept light curves 
whose shapes slightly deviate from the Paczyhski form, ei- 
ther because of a real deviation in the microlensing signal, 
as is the case for the mic rolensing event PA-99-N2 dis- 
cussed bv lAn et alJ l|2004b|) . or because the signal may be 
disturbed by artefacts or by some nearby variable stars. 

Another crucial element in the selection is the choice 
for the required sampling along the bump. In fact, while 
a good sampling is needed in order to meaningfully char- 
acterise the detected variation, demanding too much in 
this respect could lead us to exclude many bona fide can- 
didates. Using the values of ti/ 2 and to determined in the 
preceding step, we define 4 time intervals around the time 
of maximum magnification to'. [to - 3fi/2, to — f 1/2/2], [to — 
ti/a/2, to], [to, to+ti/a/2] and [to+ti/a/2, t + 3t 1/2 }. As 
a fourth cut we demand that a minimum number of ob- 
serving epochs n m i n occur in each of at least 3 of these 

2 Note that, even if it does not contain any astrophysical 
information, we must include the background pixel flux as a 
parameter in the fit to take into account its statistical fluctua- 
tion when we estimate the parameters of the Paczyriski curve. 
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Fig. 2. Top: Distribution of flux deviations at maximum 
for the selected events after the sampling cut. Bottom: 
Duration distribution for the selected events after the cut 
on R(A$). 
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Fig. 3. r and i filter 3-year light curves for 2 selected 
variations before the last cut. Upper panels, PA-99-N1, 
an accepted candidate. Lower panels: a rejected candidate. 
The dashed line is the best-fit for a Paczyhski bump with 
a sinusoidal background. The abscissae are time in days 
(JD-2451392.5). The ordinates are flux in ADU/s. 
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criterion 


number of selected light curves 


cluster detection (Q > 100) 


~ 10 s 


signal to noise ratio (Li > 40) and second bump (L2/L1 < 0.5) 


~ 4 ■ 10 4 


shape analysis : x^/dof < 10 (7 parameter Paczyriski fit) 


~ 3 ■ 10 4 


time sampling along the bump 


~ 10 4 


flux deviation: R(A$) < 21 
time width: t 1 / 2 < 25 days 


~ 1.5 ■ 10' 5 
9 


second bump analysis 


6 



Table 1. Summary of the selection criteria and number of the selected light curves. 





PA-99-N1 


PA-99-N2 


PA-00-S3 


PA-00-S4 


a (J2000) 


00h42m51.19s 


00h44m20.92s 


00h42m30.27s 


00h42m29.98s 


S (J2000) 


41°23'56.3" 


41°28'44.8" 


41°13'00.6" 


40°53'46.1" 


AO 


7' 53" 


22'04" 


4'06" 


22'33" 


ti/2 (days) 


1.83+°;- 


22.161°;- 


9 qns+ u - UV4 

iiOUJ -0 .062 


1 Q6+ u ° 9 
r.ao_ 10 


R(A&) 


20.83 ±0.10 


19.10 ±0.10 


18.80 ±0.20 


20.7 ±0.20 


V-R 


1.2 ±0.2 


1.0 ±0.1 






R-I 






0.6 ±0.1 


0.0 ±0.1 


t (JD-2451392.5) 


13.85 ±0.05 


71.70 ±0.10 


458.40 ± 0.02 


488.90 ± 0.07 


*E (days) 
u 

<p* (ADU/s) 
f g (ADU/s) 
0* (ADU/s) 
^4.max 


q+4.5 
0.0 2 7 

n7n+° ; ° 46 

U.U(U_ 030 
1 17 +0.76 
L - L ' -0.49 

0-35+ ,;? 4 
14.3±g;t 


71 1+ 4 - 1 

0.10141™ 
10.87l°;S 

q C7+O.28 
°- J '-0.25 

Q q+0-68 
a ' y -0.65 


10-4+O3 

0.070tte 
8.9+J;? 

11.71+° 

1431*;? 


0.0042+°^ 056 
0.1ll^ i5 

0.07+°^ 10 
200l?'f° 


xVdof 


1.1 


9.3 


2.1 


0.9 



Table 2. Main characteristics of the four already published microlensing candidates. AO is the projected separa- 
tion from the centre of M31. The magnitudes correspond to the maximum flux deviation and are given in standard 
Johnson/Cousins system. The results reported here are the results of the Paczyriski fit alone, even when extra infor- 
mation is available, as is the case for PA-99-N1 and PA-99-N2. 



time intervals. Clearly n m i n cannot be as large for short 
events as for long ones. We choose n min = 1, 2 and 3 for 
ti/2 < 5, t\/2 G (5,15) and t\/2 > 15 days, respectively. 
Furthermore, neither of the external intervals should fall 
at the beginning or end of one of the three seasons and at 
the same time be empty. 

The cuts described above reduce our sample of poten- 
tial events to ~ 10 4 , about one tenth of the initial set of 
selected variations, but still mostly variable stars. 

In this paper, we restrict attention to bright 
microlensing-like variations, in particular we demand 
i?(A$) < 21, although the observed deviations extend 
down to i?(A$) ~ 24 (Fig. EJ. This reduces our set of 
candidates by another factor of ~ 10. 

The Monte Carlo (Sect. 0J predicts most of the mi- 
crolensing events to be rather short. On the other hand, 
the observed ti/ 2 distribution shows a clustering of long 
variations centred on ti/ 2 ~ 60 days, most of which are 
likely to be intrinsically variable objects, and a much 
smaller set of short-duration variations (Fig. [3Jl. We de- 
mand ti/2 < 25 days, which leaves us with only 9 
Paczynski-like flux variations. 

Out of the 9 variations selected above, 5 show a sig- 
nificant second bump. We want to exclude variable stars, 



while keeping real microlensing variations that happen to 
be superimposed on a variable signal. For most variable 
stars the secondary bump should be rather similar but 
not identical to the detected one. To make use of this fact 
we perform a three-colour fit, modelling the light curve 
with a Paczyriski bump plus a sinusoidal signal, and then 
compare the time width and the flux variation of the si- 
nusoidal part with those of the Paczyriski bump. Because 
our model is very crude and because we know that vari- 
able stars may show an irregular time behaviour, we do 
not ask for a strict repetition of the bump along the base- 
line to reject a variation. We exclude a light curve if both 
the i?(A$) difference between the two bumps is smaller 
than 1 magnitude and the time width of the sinusoidal 
part is compatible with that of the bump within a fac- 
tor of 2. Three out of nine variations are excluded in this 
step. For all three the detected bump is relatively long 
(t 1/2 > 20 days) and faint (i?(A$) > 20.5). Furthermore, 
on the images the position of the second bump appears to 
be consistent with that of the detected bump, clear evi- 
dence in favour of the intrinsically variable origin of these 
variations. Two other light curves are retained, although 
they show a significant secondary bump; in both cases, 
the secondary bump is much longer than the main one. 
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Fig. 4. 3-year light curves of the four microlensing events PA-99-N1, PA-99-N2, PA-00-S3 and PA-00-S4. For each 
event, the top panel shows the whole light curve in the r filter, while the 2 lower panels display zoomed light curves 
in all bands for which data are available. Dashed lines are best-fit Paczyhski curves. The abscissae are time in days 
(JD-2451392.5). The ordinates are flux in ADU/s. 



Besides, in both cases the direct inspection on the images 
reveals that the position of the second bump is different 
from that of the detected one. In order to make clear the 
sense of the present criterion, we show (Fig. [3J the re- 
sult of the Paczyhski fit superimposed over a sinusoidal 
background for two variations. In the upper panels is an 
accepted candidate, for which the short and bright bump 
at to ~ 13 (JD-2451392.5) is clearly distinct from the un- 
derlying variable signal. In the lower panels is a rejected 
candidate. The Paczyhski signal originally selected with 
peak at to ~ 480 (JD-2451392.5) is clearly undistinguish- 
able from the underlying variable background. 



We are now left with our final selection of 6 light curves 
showing an achromatic, short-duration and bright flux 
variation compatible with a Paczyhski shape. We denote 
them PA-99-N1, PA-99-N2, PA-00-S3, PA-00-S4, PA-00- 
N6 and PA-99-S7. The letter N(S) indicates whether the 
event lies in the north (south) INT WFC field, the first 
number (99, 00, or 01) gives the year during which the 
maximum occurs, and the second has been assigned se- 
quentially, according to when the event was identified. 

In Table ^ we report in sequence each step of the 
pipeline with the number of the selected candidates re- 
maining. 
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3. Microlensing events 

3.1. POINT-AGAPE 3 years analysis results 

In this section we look at the 6 selected candidates in 
detail. In Table [5] and Figure 0] we recall the character- 
istics and light curves of the four already published can- 
didates 3 , while Table |21 and Figures [S] and [5] are devoted 
to the two new ones. The errors in i?(A$) and the colour 
index are dominated by the uncertainty in the calibration 
of the observed flux with respect to the standard magni- 
tude system, except for PA-00-N6. When the 7-parameter 
Paczyhski fit does not converge properly, the time width 
and the flux increase are estimated from a degenerate fit 
<Oouldlll99d) . 





PA-00-N6 


PA-99-S7 


a (J2000) 


00h42ml0.70s 


00h42m42.56s 


S (J2000) 


41°19'45.4" 


41°12'42.8" 


AO 


7' 16" 


3' 28" 


ti/2 (days) 
i?(A$) 


1 77 +U.57 
L - ' ' -0.60 

20.78t°;£ 


^■ xu -0.73 

20.80 ±0.10 


V -R 




0.79 ±0.14 


R-I 


U.01_ .43 




to 


491.30 ±0.07 


65.21 ±0.14 


t E (days) 
u 

4>* r (ADU/s) 
4>t (ADU/s) 
Amax 


8.3+iV 
0.07i°;« 2 
1 40+ 26 

i 7+3.2 
L - ' -1.2 
i 4 +26 




x7dof 


1.0 


1.3 



Table 3. Main characteristics of the two new microlensing 
candidates. The parameters are the same as in Tabled 



The source star of PA-99-N1 has been identified on 
HST archival images l|Auriere et all l200l[) . Fixing the 
source flux at the observed values, cf>* = 1.02 ADU/s and 
4>* g = 0.28 ADU/s, we obtain t E = 9.20 ± 0.61 days and 
uo = 0.060 ± 0.005, compatible within ltr with the values 
reported in Tabled obtained from our data alone. Finally, 
th e HST data allow u s to estimate the colour (R—I) ~ 0.9. 
In lAn et alJ (l2004bf) . we have demonstrated that PA- 
99-N2, which shows significant deviations from a simple 
Paczyhski form, is compatible with microlensing by a bi- 
nary lens. The binary-fit parameters are characterised by a 
longer time scale and higher magnification than the point- 
lens fit. In the best-fit solution we find tg = 125.0 ± 7.2 
days, m = (3.60 ± 0.37) x 10" 2 , (j>* = 4.76 ± 0.34 ADU/s, 
and a lens mass ratio ~ 1.2 x 10~ 2 . Under the assump- 
tion that the lens is associated with M31 (rather than the 
MW), the lower bounds on the angular Einstein radius 
(#e > 25 /xas) deduced from the absence of detectable 
finite-source effects implies that the source-lens relative 
velocity is v± > 280 km/s, and the source-lens distance is 



3 F u ll details can be fo und in lPaulin-Henriksson et all i200: 
I2003D : lAn et~ai1 J2004bl) . 
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Fig. 5. 3-year light curves of the microlensing event PA- 
00-N6. Panels and symbols as in Figure 
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Fig. 6. 3-year light curves of the microlensing event PA- 
99-S7. Panels and symbols as in Figure^ 



di s > 45 kpc(M/Af Q )~ 1 , where M is the lens mass. These 
facts, together with PA-99-N2's large distance from the 
M31 centre (~ 22') make it very unlikely to be due to 
an M31 star, while the prior probability that it is due to 
a MW star is extremely low. Hence, PA-99-N2 is a very 
strong MACHO candidate (either in M31 or the MW). 
The sampling and the data quality along the bump are 
also good enough to permit a reliable estimate of all 7 
parameters of the Paczyhski fit for the event PA-00-S3. 
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Fig. 7. i?(A$) — (R — I) colour-magnitude diagram for 
the ~ 10 000 variations selected before the cut on the flux 
deviation at maximum. Superimposed we show the posi- 
tions of the 6 selected candidates. The R — I colours for 
PA — 99 — N2 and PA - 99 - 57 are estimates derived 
from the observed V — R colours. 



For PA-00-S4 we obtain only a reliable lower limit on 
and accordingly an upper limit on uq, as indicated by the 
question marks in Tabic [21 

For PA-00-N6, the data allow us to evaluate the full set 
of Paczyhski parameters. Note the rather short Einstein 
time, ~ 10 days, similar to those of PA-99-N1 and PA-00- 
S3. 

As in the case of PA-99-N1 l|Paulin-Henriksson et alJ 

120031) . PA-99-S7 lies near (within 4 pixels) of a long-period 
red variable star. This induces a secondary bump, which 
is particularly visible in the i light curve. PA-99-S7 has 
been accepted by the last step of our selection pipeline, 
despite this second bump being responsible for poor sta- 
bility of the baseline. In this case, the data do not allow 
us to break the degeneracy among the Paczyhski param- 
eters and therefore do not allow a reliable estimate of the 
Einstein time. 

A colour-magnitude diagram of the ~ 10 000 varia- 
tions selected after the sampling cut is shown in Figure 
[7| Superimposed we indicate the position of the 6 varia- 
tions finally selected after all cuts. In particular, we note 
the peculiar position of PA-99-N2, which (together with 
PA-00-S3) is unusually bright relative to the other vari- 
ations. Recall that PA-99-N2 is also the longest selected 
variation, with t]/ 2 ~ 22 days. As we have already ex- 
cluded short-period variables, the sample shown is domi- 
nated by red, long-period variables of the Mira type with 
i?(A$) > 21, {R-I) > 1. For a detailed discussion of the 
variable star pop ulations detected within our dataset see 
lAn et alJ l|2004a[) . 



The spatial position for the detected events projected 
on the sky is shown, together with the INT fields, in Figure 
n Note the two new events are located within a rather 
small projected distance of M31's centre. 

3.2. Variable Contamination 



Probably the biggest single problem in the interpreta- 
tion of microlensing events drawn from faint sources is 
the possibility that the sample may be contaminated 
with rare variables. For relatively bright sources, such as 
those being detected by the thousand toward the Galactic 
bulge llUdalskil 20031. microlensing events are easily distin- 
guished from variables by their distinct shape. However, 
as the S/N declines, such discrimination becomes more 
difficult. Experiments toward the LMC provide sober- 
ing confirmation of the legitimacy of this concern. Both 
of the original microl ensing candidates rep orted by the 
EROS collaboration l|Aubourg et alJ Il993|) were subse- 
quently found to be variable stars, while some candi- 
dates found by the MACHO collaboration ( Alcock et alJ 
11 997l l2000h were also subsequently recognized as possi- 
ble or certain variab les. The SuperMACHO collaboration 
I Becker et aT1l2004|) . which probes about 2 mags fainter 
than MACHO or EROS in its microlensing search toward 
the LMC, has so far found it extremely difficult to dis- 
tinguish between genuine microlensing events and back- 
ground supernovae (C. Stubbs 2005, private communica- 
tion). Thus, when reporting a handful of microlensing can- 
didates drawn from 3 years of monitoring of a large frac- 
tion of an entire L* galaxy, we should cautiously assess 
the possibility of variable contamination. 

If variables were contaminating our sample, they would 
have to reside either in the MW or in M31 itself, or they 
could be background supernovae. We consider these loca- 
tions in turn. 

There are three arguments against MW variables: dis- 
tribution on the sky, absence of such variables in the 
Galactic microlensing studies, and lack of known classes 
of Galactic variables that could mimic microlensing. First, 
of the 5 microlensing candidates that enter our event-rate 
analysis (i.e., excluding the intergalactic microlensing can- 
didate PA-00-S4), 4 lie projected in or near the M31 bulge. 
This strongly argues that they are, in their majority, due 
to M31 sources, which are also heavily concentrated in this 
region. By contrast, Galactic variables would be spread 
uniformly over the entire field. Of course, this does not 
rule out the possibility of minor contamination by such 
variables. 

However, if there were a class of variables that could 
even weakly mimic short microlensing events with flux 
variations corresponding to i?(A<I>) < 21, then these 
would have easily shown up in Galactic microlensing ex- 
periments. For example, the OGLE-III microlensing sur- 
vey covers over 50 deg 2 toward the Galactic bulge, more 
than 100 times larger than our survey toward M31. The 
OGLE survey does not go as deep as ours because their 



POINT-AGAPE: Evidence for a MACHO contribution to Galactic Halos 



9 



telescope is smaller (1.3m) and their exposure times are 
shorter (2 min) , although these factors are somewhat com- 
pensated by their denser temporal coverage. Ignoring this 
shallower depth for the moment, and restricting consid- 
eration to < 3 kpc (where most of our foreground MW 
disc stars lie) the projected density of disc stars is about 
10 times higher in the OGLE fields than in ours because 
they lie at lower Galactic latitude. Hence, one would ex- 
pect of order 1000 times more such variables to appear 
in the OGLE fields than in ours. Of course, the major- 
ity of these would be i?(A$) ~ 21 and so of such low 
signal-to-noise ratio that they would not appear as OGLE 
candidates, or if they did, would escape recognition as 
variables. However, ~ 1/125 would lie 5 times closer and 
so be 3.5 mag brighter, i.e., i?(A$) < 17.5, correspond- 
ing to I 17, and these would have good signal-to-noise 
ratio. No such variable population is reported. A similar 
argument applies to Galactic halo stars, which would also 
be much denser in the OGLE-III fields than in ours. 

Third, there are no known candidate classes of Galactic 
variables that could mimic the M31 microlensing events. 
The one possibility is dwarf-novae, which have been re- 
ported as faking m icrolensi ng events toward the LMC 
i|Ansari et alJl995|) and M22 l)Bond et all2005|) . However, 
with typical peak absolute magnitudes of My ~ 2 
( Warnc r 1995), they would have to lie well outside the 
Galaxy to appear as i?(A$) ~ 21 fluctuations. 

While the case against M31 variables is not as airtight 
as against Galactic ones, it is still quite strong. The basic 
argument is that if the sources are in M31, then they must 
suffer luminosity changes corresponding to Mr < —3.5 on 
quite short timescales (ti/2 < 5 days for all candidates 
except PA-99-N2). There are no known classes of vari- 
ables that do this except for novae. However, novae show 
brighter variations and strongly asymmetric light curves 
characterized by slow descents (a se lection of novae var i- 
ations in our dataset is discussed in lDarnlev et alJl2004h . 
While in principle our microlensing candidates could be 
due to some new, so far unrecognized (nor even conjec- 
tured) type of stellar variability, the great brightness and 
very short timescale of the observed events impose severe 
restrictions on candidate mechanisms of variability. 

Novel mechanisms to explain the sixth event, PA-99- 
N2, would be less constrained because it is much longer, 
ti/2 ~ 22 days. However, being long as well as very bright 
(i?(A$) ~ 19), its signal-to- noise ratio is quite high. This 
permits us to check for achromaticity with very good pre- 
cision. Even the deviations from a simple Paczyhski shape 
are ac hromatic and can be reproduced by a binary-lensing 
curve l)An et al.ll2004tJ) . That is, PA-99-N2 is an excellent 
microlensing candidate on internal evidence alone. 

Finally, we remark on supernovae which, as noted 
above, plague the SuperMACHO project and also were a 
difficult contaminant for the MACHO and EROS projects. 
There are two principal arguments against supernovae. 
First, the FWHMs of all but one of the events are too 
short for supernovae while, as we have just argued, the 
sixth event is achromatic and fit by a binary-lens light 



curve and therefore almost certainly microlensing. Second 
supernovae cannot be responsible for the majority of the 
events because the supernovae would be uniformly dis- 
tributed on the sky while the actual events are highly 
clustered near the centre of M31. 

For completeness, we address one other concern related 
to variability: the possibility that the source displays a sig- 
nature of variability away from the microlensing event. In 
this case, one might worry that this "event" is actually an 
outburst from an otherwise low-level variable. Recall that 
our selection procedure actually allows for a superpixel 
to show lower-level variability in addition to the primary 
"event" that is characterized as microlensing, and to still 
be selected as a candidate. This is necessary because about 
15% of pixel light curves within 8' of the M31 centre (a re- 
gion containing most of our events) show variable-induced 
"bumps" with likelihood L\ > 40. So we would lose 15% 
of our sensitivity if we did not try to recover microlensing 
events with such secondary bumps. One event (PA-99-N1) 
out of four in this region displays such a severe secondary 
bump. This 25% rate is within Poisson uncertainties of the 
15% expectation. In addition, a second event (PA-99-S7) 
displays a secondary bump at less than half this threshold. 

It must be stressed, however, that through a Lomb 
analysis we find that neither of the source stars for these 
two events shows any sign of variability apart from the 
microlensing event. In both cases, the source of the lower- 
level variation lies several pixels from the microlensing 
event. 

In brief, while we cannot absolutely rule out non- 
microlensing sources of stellar variability, all scenarios that 
would invoke variability to explain our candidate list are 
extremely constrained, indeed contrived. 



3.3. A likely binary event 

Our selection pipeline is deliberately biased to reject flux 
variations that strongly differ from a standard Paczyhski 
light curve. In particular, it cannot detect binary lens 
events with caustic crossing. We discuss here a blue flux 
variation (R — I ~ 0) that failed to pass the x 2 cut, but 
is most probably a binary lens event: PA-00-S5. The light 
curve, which involves a short (ii/2 ~ 2 days) and bright 
peak followed by a plateau, is suggestive of binary lensing 
with a caustic crossing. The photometric follow-up of this 
event is tricky, particularly in the i band, because a faint 
resolved red object lies about 1.5 pixels away. To overcome 
this difficulty, we have used a more refined difference im- 
age photometry that includes modelling the PSF. 

We have found a binary lensing solution that convinc- 
ingly reproduces the shape of the bump. The correspond- 
ing light curve, superimposed on the data obtained us- 
ing difference image photometry, is displayed in Figure [S] 
where we show the full r light curve, zooms of the bump 
region in the r and i bands, and the ratio of flux increases 
A$ r /A$i. 
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Fig. 8. A binary solution superposed on the differential- 
photometry light curve of the binary-lens candidate PA- 
00-S5. Upper panel: full r light curve; middle panels: r and 
i zooms around the bump region, the dotted line shows 
the baseline; bottom panel: the colour ratio A$ r /A$i, 
the dash-dotted line being the average colour ratio. The 
abscissae are time in days (JD-2451392.5), the ordinates 
of the three upper panels are flux in ADU/s. 



This solution is a guess, neither optimised nor checked 
for uniqueness. The parameters are as follows: the distance 
between the two masses is d = 0.63 in unit of the Einstein 
radius Re, the mass ratio is q = 1/2; the distance of 
closest approach to the barycentre, u$ = 0.17, is reached 
at to = 411(JD-2451392.5); the Einstein time scale is ts = 
50 days; the source crosses the binary axis at an angle of 
58.5°, outside the two lenses and close to the heavy one. 

The location of PA-00-S5 is a = 00h41ml4.54s, 6 = 
40°48'37.7", J2000, some 32' away from M31's centre. This 
event cannot enter the discussion of the following sections 
because it does not survive our full selection pipeline and 
because the possibility of caustic crossings is not included 
in the simulation. Nevertheless, if this event is due to mi- 
crolensing, the lens is most probably a binary MACHO. 

3.4. Comparison with other surveys 

The first microlensing candidate reported in the direction 
of M31, AGA PE-Z1, was detected in 1995 by the AGAPE 
collaboration l|Ansari et al.lll999(l . AGAPE-Z1 is a very 
bright event, R = 17.9, of short duration, t 1 / 2 =5.3 days, 
and located in the very central region of M31, at only 
~ 42" from the centre. 

The MEGA collaboration has presented results from a 
search for microlensing events using the first 2 years of th e 
same 3-year data set analyzed here ijde Jong et alJl2004]) . 
but a different technique. In contrast to the present analy- 
sis, they do not impose any restriction on t\t-x and i?(A$). 
As a result, they select 14 microlensing candidates. All of 
them belong to our initial catalogue of flux variations. 
However, beside MEGA-7 and MEGA-11 (corresponding 
to PA-99-N2 and PA-00-S4, respectively), the remaining 
12 flux variations are fainter than allowed by our magni- 
tude cut (i?(A$) < 21). Moreover, MEGA-4, MEGA-10, 



MEGA-12 and MEGA-13 have time widths longer than 
our threshold of 25 days. 

The WeCAPP collaboration, using an original set 
of data acquired in the same period as our campaign, 
reported the detec tion of two microlensing candidates 
i Riffeser et all2003h . The candidate WeCAPP-GLl is PA- 
00-S3. We did not detect the candidate WeCAPP-GL2 
(short enough but probably too faint to be included in 
our selection) because its peak falls in a gap in our obser- 
vations. 

The NA INITAL survey has recently reported 
l),Ioshi et al.ll2005h the discovery of a microlensing candi- 
date toward M31, quite bright (i?(A$) = 20.1) but too 
long (ti/2 ~ 60 days) to be sel ected within our pipeli ne. 

Recently we have reported (jBelokurov et alJl2005|) the 
results of a search for microlensing events obtained using 
a different approach. Starting from a different catalogue of 
flux variations and using a different set of selection criteria 
(in particular, we did not include any explicit cut in t\/2 
or i?(A$)), we reported 3 microlensing candidates: PA- 
00-S3, PA-00-S4 and a third one, which is not included in 
the present selection. It is a short, bright, rather blue flux 
variation (t 1/2 = 4.1 days, i?(A$) = 19.7, R- I = 0.0), 
detected in the third year (t Q = 771 (JD-2451392.5)). In 
the present analysis it is rejected because it fails to pass 
the sampling cut: it does not have enough points on the 
rising side to safely constrain its shape. The position of 
this event, (a=00h42m02.35s, S = 40°54'34.9", J2000), 
rather far away from the centre of M31 (AO = 22'59"), is 
consistent with its being a MACHO candidate. However, 
because it does not survive the present selection pipeline, 
we do not include it in the following discussion. A further 
analysis in which w e follow a still diffe rent approach is 
currently underway ijTsanras et allEool^ . 

4. The Monte Carlo analysis 

The Monte Carlo attempts, for a given astrophysical con- 
text, to predict the number of events expected in our 
experiment, trying to mimic the observational conditions 
and the selection process. Because these can only partially 
be included in the Monte Carlo, the full simulation of our 
observation campaign must involve the detection efficiency 
analysis which is described in Sect. El 



4.1. The astrophysical model 
4.1.1. The source stars 

Source stars are drawn according to the targe t M31 lu- 
minosity profile as modelled by iKentl l|l989h . The 3- 
dimensiona l distribution of bulge stars is also taken from 
IKentl 1 19891) . The distance z of disc stars to the disc plane 
follow a 1/ c osh 2 (z/ff) distribu tion with H = 0.3 kpc as 
proposed bv lKerins et all l)200l[, . 

The colour-magnitude distributions of disc and bulge 
stars are supposed to have the characteristics of the Milky 
Way disc and bulge populations. The distribution of disc 
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stars is taken from the solar neighbourh ood data ob- 
tained by Hipparcos l|Perrvman et alJll997|) . corrected at 
the bright end for the completeness volume 4 and incor- 
porating at low lumin osity (needed for n ormalisation) a 
Besancon disc model l|Robin et al.ll2003l). For the bulge 
we again use a Besancon model llRobin et al.ll2003|) com- 
pleted at the faint end using lHan fc Gouldl l)l996(l .~We con- 
struct two distinct types of "colour-magnitude diagrams" 
(CMDs) from the Monte Carlo and show these in Figure 
El with the position of the actual detected events super- 
posed. The first is a standard CMD, which plots apparent 
magnitude versus colour for the sources of all the sim- 
ulated microlensing events that meet our selection crite- 
ria. In fact, however, while the colours and magnitudes 
of all selected-event sources are "known" in the Monte 
Carlo, they cannot always be reliably extracted from the 
actual light curves: the colours are well-determined, but 
the source magnitudes can only be derived from a well- 
constrained Paczyhski fit (while some events have only 
degenerate fits). We therefore also show in Figure El a sec- 
ond type of CMD, in which the ordinate is the magnitude 
corresponding to maximum flux increase during the event 
(i?(A$)). It is always well-determined in both the Monte 
Carlo and the data. 

To take into account the effect of the finite size of 
stars, which can be important for low mass MACHOs, 
we have to evaluate the source radii. To this end, we use 
a col our temperature re lation evaluated from the mod- 
els of iRobin et alJ l|2003[K and we evaluate the radii from 
Stefan's law u sing a table of bolometric corrections from 
iMurdiul |)2001J) . 

We did not take into account possible variations of the 
interstellar extinction across the field, although there are 
indica tions of higher extinction on the near side ijAn et alJ 
2004a). The best indicator we have of differential extinc- 
tion is the asymmetry of the surface brightness map, and 
this gives a flux attenuation by dust on the near side of 
about 10%. This is also the order of magnitude of the 
average extinction one would obtain assuming that the 
M31 disc absorption is about twice that of the MW disc. 
Indeed, as dust is confined in a thin layer, extinction only 
significantly affects the stars on the back side. Clearly an 
attenuation of about 10% would not significantly affect 
the results presented here. 
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Fig. 9. The colour/magnitude event density distribution 
predicted by the Monte Carlo. Top panel: R magnitude 
of the source star. Bottom panel: R magnitude of the flux 
increase. The observed events are superposed on the dia- 
grams. Only those events for which the source magnitude 
can be reliably extracted appear in the left panel. The 
colour scale shows the event density (in arbitrary units). 



4.1.2. The lenses 



The lenses can be stars or halo objects, with the latter 
being referred to as "MACHOs" . The stellar lenses can be 
either M31 bulge or disc stars 5 . 



4 The luminosity fu nction obtained in this w ay fully agrees 
with that presented in lJahreifi fc Wielenl jl997tl . 

5 We do not include lensing of M31 objects by stars of the 
MW disc. This can be at most of the same order of magnitude 
as M31 disc-disc lensing, which is included but turns out to be 
small. 



In the case of the bulge, we shall consider the mi- 
crolensing contribution of bulge stars with a standard stel- 
lar mass-to-light ratio. Such models form the only true lit- 
mus test for whether or not dark matter must be invoked, 
since the dark matter solution is classically required to 
explain observations which cannot be accounted for by 
known populations. The only dynamical requirement for 
our stellar bulge models is that their dynamical contribu- 
tion does not exceed the observed inner rotation curve. 
They do not need to fully reproduce the inner rotation 
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curve, though their failure to do so must be seen as evi- 
dence in itself for dark matter. We shall from here onwards 
use the term stellar bulge to denote the contribution to the 
bulge from ordinary stars. We use the term bulge by it- 
self to mean the entire dynamical bulge mass, which must 
include the stellar bulge but which may also comprise ad- 
ditional mass from unknown populations. We implicitly 
assume that the total bulge mass is fixed by the rotation 
curve. We set out here to discover whether or not the rate 
predicted by known stellar bulge and disc populations can 
feasibly account for our observed microlensing candidates. 

Disk stellar lenses The d isc mass distribution is the same 
as in lKerins etafl (|200l|) : 



p exp(-J)/cosh 2 (|) 



with po = 0.3M Q pc~ 3 , H = 0.3 kpc and h — 6.4 kpc. 

The mass of the disc is 3 10 10 M©, corresponding to an 
average disc mass-to- light ratio M/L B about 4. 

Bulge stellar lenses The bulge 3-dimensional mass distri- 
bution is taken to be proportional to the 3-dimensional lu- 
minosity distribution, which means that the bulge (M/L) 
ratio is position independent. Assuming that the M31 stel- 
lar bulge is s imilar to that o f the Milky Way, one can es- 
timate from lHan fc Gouldl l)2003l) that M/L B ~ 3 and 
that it cannot exceed 4 (corresponding to bulge masses 
of 1.5 and 2 10 10 M within 4 k pc). This can also b e in- 
ferred by combining results from|Zoccali et al.l |2000) and 
iRoeer et al.l l|l986|) . lHan fc Gouldl (|2003[) have shown that 
this stellar M/L accurately predicts the optical depth that 
is observed toward the MW bulge. 

Estimates higher than the above values for the 
total bulge and dis c M/ L b ha v e been quot ed on 



dynamical grounds llKentl Il989t iKerins et alJ 
" Baltz et al.ll2003UWidrow et al.ll2003HGeehan etafl 



2001 



2005 



Widrow fc Dubinskil |200 5 | ) and used t o make predictions 
on self lensing (e.g. iBaltz et al.ll2003|) . In these dynami- 
cal studies a heavy bulge (M ~ 410 10 M©, M/L B ~ 8) 
is typically associated with a light disc (M ~ 3 10 10 Mq, 
M/Lb ~ 4), whereas a light bulge (M ~ 1.510 10 M Q , 
M/L B ~ 3) goes with a heavy disc (M ~ 710 10 M Q , 
M/Lb ~ 9). As stated above, such large M/Lb ratios 
mean that some kind of dark matter must be present as 
no known ordinary stellar populations can provide such 
high M/Lb ratios. We shall refer to these solutions to 
evaluate upper bounds on the self-lensing contribution in 
Sect. El 

The stellar mass function is taken from iKer ins et alJ 



(2001 

dN 
dm 



,-0.75 



(0.08 M Q < m< 0. 5M ) 

(0.5 Mq < TO < 10 Mq) 



(4) 



The corresponding average stellar mass is < m >~ 0.65 
M©. We ha ve also considered st eeper mass functions, as 
proposed bv lZoccali et al.l l|2000|) . for which < to >^ 0.55 



M , or bv lHan fc Gouldl J2003I) . for which < m >~ 0.41 
M . Our results turn out to be rather insensitive to this 
choice. 

Halo lenses (MACHOs) The MW and M31 halos are mod- 
elled as spherical nearly isothermal distributions with a 
core of radius a : 



P\ r ) — ~5 o 

' K ' a 2 + r 2 



(5) 



The central halo density is fixed, given the core radius, to 
produce the asymptotic disc rotation velocity far from the 
galactic centre. For the Milky Way the core radius clmw 
is chosen to be 5 kpc. For M31 we choose aM3i = 3 kpc for 
our reference model but we have also tried olmh = 5 kpc. 
A larger value for the core radius decreases the number 
of expected events and makes their spatial distribution 
slightly less centrally concentrated. 

As nothing is known about the mass function 
of putative MACHOs, we try a set of single val- 
ues for their masses, ranging from 10~ 5 to 1 Mq 
(10~ 5 , 10~ 4 , 10- 3 , 10~ 2 , HT 1 , 0.5 and 1 M Q ). We shall 
refer to these as "test masses" . 

4.1.3. Bulge geometry 

The most important contribution to self lensing comes 
from stellar bulge lenses and/or stars. As the event 
rates are proportional to the square root of the lens- 
source distan ce, th e bulge geometry may play an impor- 
tant role. In iKentl Jh)M>). the bulge is described as an 
oblate axisymmetric ellipsoid, and the luminosity den- 
sity is given as a function of the elliptical radius r e — 
y ' x 2 + y 2 + (z/(l — e(r e ))) 2 , where z is the distance to 
the M31 plane and e(r e ) is the ellipticity, which varies as 
a function of the elliptical radius, r e . The Kent bulge is 
quite flattened, and one may wonder if a less flattened 
model would result in more self-lensing events. To check 
this, we have run the Monte Carlo for a spherical bulge 
(e = 0), keeping the total bulge mass and luminosity fixed. 
The expected number of both bulge-disc and disc-bulge 
events rise both by about 10%. On the other hand, in 
absolute terms, the more numerous contribution of bulge- 
bulge events decreases by about 5% for a net total increase 
of ~ 2%. That is, the substitution of a spherical bulge for 
an elliptic one has almost no impact on the total rate of 
stellar bulge lensing. This can be traced to the fact that 
M31 is seen nearly edge on, which reduces the impact of 
distances perpendicular to the disk. 

4.1.4. Velocities of lenses and sources 

The relative velocities of lenses and sources strongly in- 
fluence the rate of microlensing events. The choice of 
the velocities adopted in our r eference mode l , here after 
called model 1 , is ins pired from lWidrow et all (|2P03) and 
lOeehan et al.1 l|2005l) . We stress that the bulge velocity 



POINT-AGAPE: Evidence for a MACHO contribution to Galactic Halos 



13 



dispersion is sensitive not to the mass of the stellar bulge 
component which contributes to the self-lcnsing rate, but 
to the mass of the entire bulge, which may additionally 
include unknown lensing populations. We have tested the 
effect of changing the bulge velocity dispersion and the 
M31 disc rotation velocity in models 2 to 5. The velocities 
of the various M31 components adopted for each model 
are displayed in Table 0] The solar rotation velocity is al- 
ways taken to be 220 km/s and halo dispersion velocities 
are always l/y/2 times the disc rotation velocities. All ve- 
locity dispersions are assumed isotropic, with the values 
given being 1-dimensional. 

To get an insight into the model dependence of the 
Monte Carlo predictions, it is useful to split the observed 
spatial region into an "inner" region where most self- 
lensing events are expected, and an "outer" region which 
will be dominated by MACHOs if they are present. We 
set the boundary between the two regions at an angular 
distance of 8' from the centre of M31. 

The effect of changing the velocities for the models 
displayed in Tabled is shown in Table [SJ This gives the 
relative change with respect to our reference model (for a 
MACHO mass of 0.5 M and ajwai = 3kpc). 

Beside these normalisation changes, the distributions 
of the number of events, as a function of ii/2, the angular 
distance to the centre of M31, and the maximum flux in- 
crease, all turn out to be almost independent of the model. 



Model 


bulge velocity 


disc rotation 


dispersion (km/s) 


velocity (km/s) 


1 (reference) 


120 


250 


2 


120 


270 


3 


120 


230 


4 


140 


250 


5 


100 


250 



Table 4. Velocities of M31 components (km/s). The bulge 
rotation velocity and disc velocity dispersion are fixed at 
40 km/s and 60 km/s, respectively. 





Self Lensing 


MACHOs 


Model 


Inner 


Outer 


Inner 


Outer 


region 


region 


region 


region 


2 


0.97 


0.98 


1.15 


1.21 


3 


0.97 


0.96 


0.84 


0.81 


4 


1.03 


1.03 


0.98 


1.01 


5 


0.92 


0.90 


0.98 


0.99 



Table 5. The velocity dependence of the number of ex- 
pected events. The numbers are the ratio of the number 
of expected events for models of Table 0] to the same 
number in the reference model (with M=0.5 M© and 
0-M31 — 3kpc). The number of events expected in the 
reference model, corrected for detection efficiency, are dis- 
played in Table El of Sect.© 



4.1.5. Consistency check 

To check the consistency of our Monte Carlo, we have 
computed the optical depths of the halo both analytically 
and with the Monte Carlo. The re sults are identical and 
consistent with p ublished results llCvuk fc Crottsl 1200(1 
iBaltz fc Silkll200(t . 

4.2. Modelling the observations and the analysis 

The Monte Carlo generates and selects light curves in- 
cluding part of the real observational conditions and of 
the selection algorithm. 

Reproducing the photometry conditions in the Monte 
Carlo is an important issue, so we use the same filter as in 
the real experiment. This is also true for the colour equa- 
tions, which relate fluxes to standard magnitudes in the 
reference image. In generating the light curves, all pho- 
tometric coefficients relating the observing conditions of 
the current image to those of the reference are used in 
the Monte Carlo, except for those related to the seeing 
correction. 

The observation epochs and exposure times reproduce 
the real ones, with one composite image per night. In or- 
der to avoid counting the noise twice, no noise has been 
added to the Monte Carlo light curves; it only enters via 
the error bars. As we further discuss in Sect. El an impor- 
tant condition for the efficiency correction to be reliable 
is that the Monte Carlo should not reject events that the 
real analysis would have accepted. For this reason, the er- 
ror bars in the Monte Carlo light curves only include the 
photon noise, and, for an event to be considered detected, 
we demand only the minimum condition that the corre- 
sponding bump rise above the noise (that is, L > 0, where 
L is the estimator introduced in equation 0. 

4.3. Event properties 

The main observational properties of the events are the R 
magnitude corresponding to their flux increase (i?(A$)) 
and their duration, which we characterise by the full- 
width-at-half-maximum of the bump, ii/2- The CMDs are 
displayed in Figure El We show in Figure UHl the expected 
distribution of i?(A$), the R magnitude of the sources and 
the expected t 1 / 2 distribution for two MACHO masses. 
The distribution of t 1 / 2 , quite concentrated toward short 
durations, has motivated our choice for the low-duration 
cutoff in the selection. 

5. Detection efficiency 

5.1. The event simulation 

The Monte Carlo described in the previous section does 
not take into account all the effects we face in the real 
data analysis. Therefore, its results, in particular the pre- 
diction on the expected number of events, can only be 
looked upon as an upper limit. In order to make a mean- 
ingful comparison with the 6 detected events, we must sift 
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t m (days), m = 0.1 M s||n 

Fig. 10. The expected distribution of the R magnitude of 
the flux increase and the source stars (upper panels). The 
expected t 1 / 2 distribution for MACHO masses of 1 M Q 
and 0.1 M Q (lower panels). 

the Monte Carlo results through an additional filter. This 
is the "detection efficiency" analysis described in this sec- 
tion, wherein we insert the microlensing events predicted 
by the Monte Carlo into the stream of images that con- 
stitute our actual data set 6 . This allows us to calculate 
the detection efficiency relative to the Monte Carlo and to 
obtain a correct estimate of the characteristics and total 
number of the expected events. 

The main weakness of the Monte Carlo in reproduc- 
ing the real observations and analysis stems from the fact 
that it only generates microlensing light curves, so that 
it cannot take into account any aspect related to image 
analysis. 

The Monte Carlo does not model the background of 
variable stars, which both gives rise to high flux varia- 
tions that can mimic (and disturb the detection of) real 
microlensing events, and generates, from the superposi- 
tion of many small-amplitude variables, a non-Gaussian 
noise that is very difficult to model. 

6 We refer to this analysis as "event simulation", not to be 
confused with the Monte Carlo simulation described in the 
previous section. 



As regards the selection pipeline itself, the Monte 
Carlo cannot reproduce the first, essential, cluster detec- 
tion step described in Sect. 12.21 Therefore, it cannot test 
to what extent the presence in the images of variations 
due to the background of variable stars, seeing variations, 
and noise, affect the efficiency of cluster detection, locali- 
sation, and separation. 

The Monte Carlo includes neither the seeing variations 
nor their correction nor the residuals of the seeing stabil- 
isation, which also give rise to a non-Gaussian noise. 

In principle, it would be possible to reproduce, within 
the Monte Carlo, the full shape analysis along the light 
curve followed in our pipeline. However, the results on the 
real data turn out quite different, mainly because the real 
noise cannot be correctly modelled analytically. 

In practice, no noise is included in the Monte Carlo 
light curves, because the full noise is already present in 
the images. Moreover, we have to be careful not to exclude 
within the Monte Carlo variations that the real pipeline 
is able to detect. As a consequence, the "shape analysis" 
in the Monte Carlo is quite basic. We demand only that 
the (noiseless) variations reach 3 a above the baseline for 
three consecutive epochs, where a includes only the pho- 
ton noise. 

The time sampling of our data set is fully reproduced 
by the Monte Carlo. However, the sampling criterion along 
the bump is only implemented in a very basic way by 
demanding that the time of maximum magnification lie 
within one of the 3 seasons observation. 

A typical Monte Carlo output involves ~ 20 000 events 
per CCD. However, adding 20 000 events per CCD would 
significantly alter the overall statistical properties of the 
original images (and therefore of the light curves). In or- 
der that the event simulation provide meaningful results, 
we cannot add that many events. On the other hand, the 
more events we add, the larger the statistical precision. 
Particular care has to be taken to avoid as much as pos- 
sible simulating two events so near each other that their 
mutual interaction hinders their detectability. Of course, 
these difficulties are worse around the centre of the galaxy, 
where the spatial distribution of the events is strongly 
peaked. Balancing these considerations, we choose to sim- 
ulate 5 000 events per CCD. The results thus obtained are 
compatible, with much smaller errors, with those we ob- 
tain by adding only 1 000 events (in which case the crowd- 
ing problems mentioned above are negligible). 

Each event generated by the Monte Carlo is endowed 
with a "weight" 7 , it;,, so when we refer to simulated events, 
"number" always means "weighted number" . Thus n s i m = 
J2i w h with statistical error An S i m — wf, where the 

sum runs over the full set of simulated events. 



7 As often in Monte Carlo simulations, a weight is ascribed 
to each generated event. This weight carries part of the infor- 
mation on the probability for the event to occur, before and 
independently of any selection in either the Monte Carlo or the 
event simulation. 



POINT-AGAPE: Evidence for a MACHO contribution to Galactic Halos 



15 



Let rib = n s + n>r be the number of events we simu- 
late on the images, where n s and n r are respectively the 
number of selected and rejected events at the end of the 
analysis pipeline. We define the detection efficiency as 

_ ns_ 
~ rib 

and the relative statistical error is then 

f Ae\ 2 (n r An s ) 2 + (n s An r ) 2 
V e J (n&n s ) 

Once we know e, we can determine the actual number of 
expected events, n exp — e , where is the number 
expected from the Monte Carlo alone. 

The event simulation is performed on the images af- 
ter debiasing and flatfielding, but before any other reduc- 
tion step. We use the package DAOPHOT within IRAF. 
First, starting from a sample of ~ 200 resolved stars per 
CCD, for each image we evaluate the PSF and the rela- 
tive photometry with respect to the reference image. Then 
we produce a list of microlensing events, randomly chosen 
among those selected within the Monte Carlo. For each 
event, using all the light curve parameters provided by 
the Monte Carlo as input, we add to each image the flux 
of the magnified star at its position, convolved with the 
PSF of the image (taking due account of the required ge- 
ometrical and photometric calibration with respect to the 
reference image). We then proceed as in the real analysis. 
In particular, after image recalibration, we run the selec- 
tion pipeline described in Sect. 12.21 In short, the scope of 
the event simulation is to evaluate how many "real" mi- 
crolensing events are going to be rejected by our selection 
pipeline. We test the event simulation procedure by com- 
paring the mean photometric dispersion in the light curves 
of observed resolved stars to those of simulated, stable, 
stars of comparable magnitude. We find good agreement. 

In the selection pipeline, it is essential to use data 
taken in at least two passbands in order to reject vari- 
able objects. Indeed, we test achromaticity with a simul- 
taneous fit in two passbands and, in the last step of the 
selection, we test whether a secondary bump is compatible 
with being the second bump of a variable signal. Here, us- 
ing i band data is important because the main background 
arises from long-period, red variable stars. 

In the event simulation, we want to evaluate what frac- 
tion of the Monte Carlo microlensing events survive the 
selection pipeline. For these genuine microlensing events, 
we expect the use of two passbands to be less impor- 
tant. In fact, microlensing events are expected to pass the 
achromaticity test easily. Moreover, because the events we 
simulate are short and bright, the microlensing bump is 
in general quite different from any possible, very often 
long 8 , secondary bump, and most simulated events pass 
the secondary-bump test. Indeed, we have checked on one 

8 Short-period variable objects have already been removed 
since they are easily recognised from their multiple variations 
within the data stream. 



CCD that we get the same result for the detection effi- 
ciency whether we use data in both r and i bands or in r 
alone. For this reason, we have carried out the rest of the 
event simulation with r-band data only. 

5.2. The results 

For each CCD (with 4 CCDs per field) we simulate at 
most 5000 microlensing events, randomly chosen among 
those selected within the Monte Carlo, and subject to 
conditions reflecting the selection cuts. We only simulate 
events that are both bright (i?(A$) < 21.2) and short 
(ti/2 < 27 days). These thresholds are looser than those 
used in the selection (i?(A$) < 21.0 and t 1 / 2 < 25 days) 
because we want to include all events that can in princi- 
ple be detected by the pipeline. These enlarged cuts reflect 
the dispersion of the difference between the input and out- 
put event parameters of the event simulation. To test this 
choice, we have also run some test jobs using slightly dif- 
ferent input cuts. For instance, if one uses the looser cuts 
i?(A$) < 21.5 and tx/2 < 30 days, the number of events 
predicted by the Monte Carlo is larger, but the efficiency 
turns out to be smaller. The two effects cancel, and the 
end result for the number of expected events corrected for 
detection efficiency remains unchanged. For each CCD we 
run the event simulation for our test masses. As in the 
real analysis, we mask the very central region of M31. 

The detection efficiency depends mainly on the dis- 
tance from the centre of M31, the time width, and the 
maximum flux increase. We run the event simulation only 
for model 1 (Sect. |4~Q|) and a M31 core radius a A/31 = 
5kpc. In fact, there is no reason for the efficiency at a 
given position in the field to depend on the core radius. 
It could in principle depend on distributions of the time 
width and the maximum flux increase, but we have seen 
that these distributions are almost model-independent. 

Finite-source effects can produce significant deviations 
from a simple Paczyhski shape, and this can be quite im- 
portant toward M31, where most sources are giant stars. 
We expect this effect to be particularly relevant for low 
mass MACHOs. The events generated by the Monte Carlo 
(Sect. 0| and entered in the event simulation include finite- 
source effects, although the microlensing fit in the se- 
lection pipeline uses only simple Paczynski curves. This 
causes an efficiency loss, which we evaluate as follows: 
we run an event simulation, for one CCD and all test 
masses, without including finite-source effects in the in- 
put events, and then evaluate the associated efficiency rise. 
This ought to be of the same order as the efficiency loss in 
the real pipeline. For masses down to 10~ 2 Mq the change 
turns out to be negligible. For masses smaller or equal to 
10~ 3 Mq, it is of the order of 20% or less. 

The detection efficiency depends on position in the 
field primarily through the distance to the centre of M31. 
At a given distance we find no significant difference be- 
tween the various CCDs. At angular distances larger than 
8' the efficiency is practically constant. In the region in- 
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criterion 


e (AO < 4') 


e (4 < AO < 8') 


e(A© > 8') 


cluster detection (Q > 100) 


46.3 ±4.1 


62.7 ± 1.5 


76.4 ± 0.4 


Li > 40 and L2/U < 0.5 
X 2 /dof < 10 
sampling 
t 1/2 < 25 days, i?(A$) < 21 
variable analysis 


40.0 ±4.0 
35.7 ±3.8 

17.1 ± 2.9 
14.7 ±2.8 
14.7 ±2.8 


57.9 ± 1.5 
54.0 ± 1.5 
31.9 ± 1.4 
25.2 ± 1.3 
25.2 ± 1.3 


72.5 ±0.4 
66.7 ±0.4 
33.7 ±0.4 
28.5 ±0.4 
28.5 ±0.4 



Table 6. Detection efficiency relative to the Monte Carlo (in percent), for a MACHO mass M= 0.5 M Q , evaluated at 
each step of the selection pipeline in different ranges of distance from the centre of M31. 



MACHO mass (M ) 


e (AO < 4') 


£(4 < AO < 8') 


£(A© > 8') 


1 


19.0 ±3.0 


24.2 ± 1.3 


29.7 ±0.4 


5 • 10 _1 


14.7 ±2.8 


25.2 ± 1.3 


28.5 ±0.4 


10- 1 


18.8 ±3.4 


22.1 ± 1.3 


26.4 ±0.4 


10 -2 


17.0 ±3.7 


21.8 ± 1.6 


23.5 ±0.5 


10 -3 


10.1 ±3.2 


14.1 ± 1.6 


15.6 ±0.5 


10- 4 


2.4 ± 1.5 


8.9 ±2.5 


9.5 ±0.5 


10 -5 


0.37 ±0.43 


5.4 ±2.2 


6.2 ±0.7 


self lensing 


17.8 ± 1.2 


22.6 ± 0.6 


26.9 ±0.3 



Table 7. Detection efficiency relative to the Monte Carlo (in percent), for our test set of MACHO masses and for self 
lensing, for the same distance ranges as in Table El 



side 8', the efficiency steadily decreases toward the centre. 
This can be traced to the increase of both the crowding 
and the surface brightness. Indeed, the drop of efficiency 
in the central region mainly comes from the first step of 
the selection pipeline, namely the cluster detection. 

Table shows the contribution of the successive steps 
of the analysis to the total loss of efficiency. The distance 
to the centre of M31 is divided into 3 ranges (AO < 
4', 4' < AO < 8' and AO > 8'). The MACHO mass is 
0.5 Mq but the qualitative features discussed below are 
the same for all masses. We have isolated the first step of 
the analysis, the cluster detection, which is implemented 
on the images, while the others are performed on the light 
curves. As emphasised earlier, the increase in crowding 
and surface brightness near the centre causes a significant 
drop of efficiency in the two central regions. Most of the 
dependence of the efficiency on the distance to the cen- 
tre arises from this step, whereas the effects of all other 
steps, acting on light curves, are nearly position indepen- 
dent. Note the loss of efficiency by almost a factor of 2 
associated with the sampling cut. This is not surprising 
as this cut is implemented in the Monte Carlo in only a 
very basic way. 

Table [7| gives the detection efficiency for our test set 
of MACHO masses after the full event selection. Down to 
a mass of 10 -2 M Q , we find no significant differences be- 
tween self-lensing and MACHO events. This reflects the 
fact that their main characteristics do not differ signifi- 
cantly on average. For very small masses, we find a drop 
in the efficiency, due to both the smaller time widths of 
the bump and finite-source effects. 



6. Results and halo fraction constraints 

In this section, we present the result of the complete simu- 
lation, the Monte Carlo followed by the event simulation, 
and discuss what we can infer about the fraction / of 
MACHOs present in the halos of M31 and the MW from 
the comparison with the data presented in Sect. [31 

In Table |H1 we present the expected numbers of self- 
lensing and halo events (for a full halo and two different 
values of the core radius) predicted by the full simula- 
tion in the three distance ranges AO < 4', 4 < AO < 8' 
and AO > 8'. The self-lensing results, given for a stellar 
bulge M/Lb ratio equal to 3, are dominated by stellar 
bulge lenses and therefore scale with this ratio. This must 
be compared with the 5 microlensing events reported in 
Sect. [3J PA-00-S4, which is located near the line of sight 
toward the M32 gal axy, is likely an intergal actic microlens- 
ing event l|Paulin-Henriksson et alJ 12002^1 and therefore 
not included in the present discussion. Accordingly, we 
have excluded from the analysis a 4' radius circular region 
centred on M32. 

The main issue we have to face is distinguishing self- 
lensing events from halo events. This is particularly impor- 
tant as the number of expected MACHO and self-lensing 
events is of about the same order of magnitude if the halo 
fraction is of order 20% or less as in the direction of the 
Magellanic clouds. 

Although the observed characteristics of the light 
curves do not allow one to disentangle the two classes 
of events, the spatial distribution of the detected events 
(Fig-Q can gi ye us useful insights. While most self-lensing 
events are expected in the central region, halo events 
should be more evenly distributed out to larger radii. In 
Figure together with the distance dependence of the 
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AO < 4' 


4 < AO < 8' 


AO > 8' 


mass (M Q ) 




Halo, OM3i = 3 kpc 
1 

5 • ltr 1 

HT 1 
1(T 2 
1(T 3 
10" 4 
10" 5 


0.70 ± 0.12 
0.81 ±0.17 
1.63 ±0.32 
1.93 ±0.45 
0.72 ±0.27 
0.064 ± 0.042 
0.002 ± 0.002 


1.38 ± 0.08 
1.93 ± 0.11 
3.00 ± 0.20 
3.85 ± 0.30 
2.20 ±0.30 
0.60 ±0.18 
0.034 ±0.015 


2.96 ± 0.04 
4.18 ±0.08 
8.10 ±0.23 
12.65 ±0.29 
9.17 ±0.29 
3.09 ±0.18 
0.42 ± 0.06 


Halo, OM3i = 5 kpc 
1 

5 • KT 1 
HT 1 
1(T 2 
10~ 3 
10" 4 
10" 5 


0.60 ± 0.10 
0.74 ±0.18 
1.30 ±0.25 
1.41 ±0.34 
0.81 ±0.30 
0.15 ±0.15 
0.002 ± 0.002 


1.11 ± 0.06 
1.57 ±0.09 
2.52 ±0.16 
3.63 ± 0.29 
2.07 ±0.26 
0.49 ±0.15 
0.048 ± 0.022 


2.48 ± 0.04 
3.63 ± 0.09 
6.94 ±0.12 
11.29 ±0.24 
8.41 ± 0.26 
2.83 ±0.16 
0.40 ± 0.05 


self lensing 


0.29 ±0.02 


0.29 ± 0.01 


0.16 ±0.01 



Table 8. The expected number of MACHO and of the self-lensing events, corrected for efficiency, for the models with 
aM3i = 3 kpc and OM3i = 5 kpc, in three different ranges of distance from the M31 centre. The stellar bulge (disc) 
M/Lb ratio is equal to 3 (4). 



fflM3i = 3 kpc om3i = 5 kpc 



Mass (M©) 


/iNF 


/max 


/sup 


/iNF 


/max 


/sup 


1 


0.27 


0.81 


0.97 


0.29 


0.97 


0.97 


5 • 10" 1 


0.22 


0.57 


0.94 


0.24 


0.67 


0.96 


HT 1 


0.13 


0.31 


0.74 


0.15 


0.37 


0.83 


10~ 2 


0.08 


0.21 


0.51 


0.09 


0.23 


0.57 


io- 3 


0.11 


0.29 


0.73 


0.12 


0.31 


0.76 


10~ 4 


0.20 


0.77 


0.96 


0.18 


0.81 


0.96 


10~ 5 


0.12 


1.00 


0.97 


0.10 


1.00 


0.97 



Table 9. Results for the halo fraction /: the 95% CL lower bound (/inf) and upper bound (/sup), an d maximum 
probability (/max) are displayed for flM3i = 3 kpc and a,M3i = 5 kpc. In both cases, the stellar bulge (disc) M/Lb 
ratio is 3 (4). 



detection efficiency, we show the expected spatial distri- 
bution of self lensing and 0.5 M© MACHO events (full 
halo) . The observed events are clustered in the central re- 
gion with the significant exception of PA-99-N2, which is 
located in a region where the self-lensing contamination 
to MACHOs events is expected to be small. 

The key aspect of our analysis is the comparison of 
the expected spatial distribution of the events with that 
of the observed ones. In order to carry out this comparison 
as precisely as possible, we divide the observed field into 
a large number of bins, equally spaced in distance from 
M31's centre. We present here an analysis with 20 bins 
of 2' width, but we have checked that the results do not 
change significantly if we use either 40 bins of 1' width or 
10 bins of 4' width. 

6.1. The halo fraction 

The first striking feature in the comparison between pre- 
dictions and data is that we observe far more events than 



predicted for self lensing alone. Therefore, it is tempting 
to conclude that the events in excess with respect to the 
prediction should be considered as MACHOs. This state- 
ment can be made more quantitative: given a MACHO 
halo fraction, /, we can compute the probability of getting 
the observed number of events and, by Bayesian inversion, 
obtain the probability distribution of the halo fraction. 

As already outlined, we bin the observed space into 
Nbin equally spaced annuli and then, given the model pre- 
dictions Xi (i — 1 ... N{,in), obtain the combined proba- 
bility of observing in each bin rii events. The combined 
probability is the product of the individual probabilities 
of independent variates nf. 



N bi n 1 

P{ni\xi) = Y[ — exp(-Xi)x?- (6) 
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Fig. 12. Most probable value, upper and lower 95% CL 
limit for the halo fraction as a function of the MACHO 
mass for aM3i = 3kpc and stellar bulge (disc) M/Lb = 
3(4). 



Fig. 11. Predictions of the full simulation as a function 
of the distance from the centre of M31. Upper pannel: 
efficiency correction (for MACHOs); central panel: ex- 
pected number of MACHO events (full halo, M = 0.5 M , 
QM31 = 3kpc); bottom panel: expected number of sclf- 
lensing events (for a stellar bulge (disc) M/Lb = 3(4)). 
The vertical lines indicate the position of the observed 
events, the dashed line corresponds to PA-00-S4, which 
has been excluded from the analysis because it is a prob- 
ably M31/M32 intergalactic event. 



For a given a model, the different not indepen- 

dent: they all depend on the halo fraction / via the equa- 
tions 



6.2. Self-lensing background ? 

The fact that 4 out of the 5 observed events lie within 8' 
from the centre of M31 could be suggestive of self-lensing 
origin, implying that we underestimate this contribution. 
However, in the Monte Carlo section we have already seen 
that the velocity dependence of our results is very weak. 
For models 2 (3), where the change is maximum, the 95% 
CL lower limit on / in the mass range 0.1-1 M Q is shifted 
by about - (+) 0.02. Furthermore, M/Lb ratios larger 
than 4 cannot be accommodated by known stellar popu- 
lations. Still, for comparison, we have considered models 
for which, on dynamical grounds, the M/Lb ratio of either 
the disc or the bulge take values up to ~ 8 — 9. One can see 
from Table 1101 that our conclusions are not qualitatively 
altered. This can be partly attributed to the occurence of 
PA-99-N2 22' away from the M31 centre. 



Xi = hif + Si, (7) 

where hi and are the numbers of events predicted in 
bin i for a full MACHO halo and self lensing, respectively. 
A model specifies h and s, so the probability depends on 
only one parameter, /. It is therefore possible to evaluate 
lower and upper limits at a given confidence level for the 
halo fraction /. 

In Figure ^] and Table |5| we display the 95% confi- 
dence level (CL) limits obtained in this configuration for 
OM3i = 3kpc and M/Lb = 3. We get a significant lower 
limit, /inf > 20%, in the mass range from 0.5 to 1 Mq. No 
interesting upper bound on / is obtained except around 
a mass of 10~ 2 M Q (/sup = 50%). We also show in Table 
El the same limits for aM3i = 5kpc. As the predicted halo 
contribution is smaller, the inferred lower limit on / is 
slightly larger. 



Bulge M/L B 


Disc M/L B 


nsL 


P(f = 0) 


/inf 
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0.72 


1CT 4 


0.22 
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1.1 


10" 3 


0.17 
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1.5 


410" 3 


0.15 



Table 10. For different sets of values of stellar bulge and 
disc M/Lb (Sect. ETl~2*|) we report the number of expected 
self-lensing events, corrected for the efficiency, the proba- 
bility for the signal to be a Poisson fluctuation for a / = 
halo and, for a M = 0.5 M MACHO population MW and 
M31 halos with ajwai = 3kpc, the 95% CL lower bound 
for the halo fraction /. 



One can also question the bulge geometry. However, we 
have seen that assuming a spherical bulge with the same 
mass and luminosity does not alter the results. One could 
also think of a bar-like bulge. This possibility has been 
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considered bv lGerhardl l|l986|) . who has shown that unless 
a would be bar points toward us within 10°, its ellipticity 
does not exceed 0.3. This cannot produce a significant 
increase of the self-lensing prediction. Even if a bar-like 
bulge points toward us and is highly prolate, it cannot 
explain event PA-99-N2. 

Clearly, unless we grossly misunderstand the bulge of 
M31, our events cannot be explained by self lensing alone. 

Still, in view of our low statistics, we could be facing 
a Poisson fluctuation. However, this is highly improbable: 
given the prediction of our simulation, the probability of 
observing 5 self-lensing events with the observed spatial 
distribution is P(f = 0) ~ 1CT 4 for a M/L B = 3 (4) M31 
stellar bulge (disc), and remains well below ~ 10~ 2 even 
for much heavier configurations ( Table ITU)l . 

7. Conclusions 

In this paper, we present first constraints on the halo frac- 
tion, /, in the form of MACHOs in the combined halos of 
M31 and MW, based on a three-year search for gravita- 
tional microlensing in the direction of M31. 

Our selection pipeline, restricted to bright, short- 
duration variations, leads us to the detection of 6 can- 
didate microlensing events. However, one of these is likely 
to be a M31-M32 intergalactic self-lensing event, so we do 
not include it when assessing the halo fraction /. 

We have thoroughly discussed the issue of the possi- 
ble contamination of this sample by background variable 
stars. Indeed, we are not aware of any class of variable 
stars able to reproduce such light curves, therefore we have 
assumed that all our candidates are genuine microlensing 
events. 

To be able to draw physical conclusions from this re- 
sult, we have constructed a full simulation of the expected 
results, which involves a Monte Carlo simulation com- 
pleted by an event simulation to account for aspects of 
the observation and the selection pipeline not included in 
the Monte Carlo. 

The full simulation predicts that M31 self lensing alone 
should give us less than 1 event, whereas we observe 5, 
one of which is located 22' away from the M31's centre, 
where the expected self-lensing signal is negligible. As the 
probability that we are facing a mere Poisson fluctuation 
from the self-lensing prediction is very small (~ 0.01%), 
we consider these results as evidence for the detection 
of MACHOs in the direction of M31. In particular, for 
«M3i = 3kpc and a M/Lg ratio for the disc and stel- 
lar bulge smaller than 4, we get a 95% CL lower limit of 
20 — 25% for /, if the average mass of MACHOs lies in 
the range 0.5-1 Mq. Our signal is compatible with the one 
detected in the directio n of the Magellanic clouds by the 
MACHO collaboration l|Alcock et al.ll2000h . 

We have also considered models that, on dynamical 
grounds, involve higher disc or stellar bulge M/ Lb ratios. 
However, because of the spatial distribution of the ob- 
served events, the conclusion would not be qualitatively 



different. Indeed, because of the presence of the event PA- 
99-N2 22' away from the M31 centre where self lensing 
is negligible, the lower bound on / would not pass below 
~ 15% even in the most extreme models considered. 

Finally, the observed events can hardly be blamed on 
the geometry of the bulge. Indeed, the number of predicted 
self-lensing events cannot be significantly increased unless 
it has a highly prolate bar-like structure exactly pointing 
toward us. However, even this improbable configuration 
would not explain one of the events, which definitely oc- 
curs outside the bulge. 

Beside the 5 events selected by our pipeline, we have 
found a very likely candidate for a binary lensing event 
with caustic crossing. This event occurs ~ 32' away 
from M31's centre, where one can safely ignore self lens- 
ing. Therefore, although included in neither our selection 
pipeline nor our discussion on the halo fraction, this de- 
tection strengthens our conclusion that we are detecting 
a MACHO signal in the direction of M31. 

To get more stringent constraints on the modelling of 
M31, better statistics are badly needed. To achieve this 
goal using our data, we plan to extend the present analy- 
sis in a forthcoming work by looking for fainter variations. 
Another option would be to lift the duration cut. However, 
we consider this less attractive, because the contamination 
by the background of variable stars would be much larger 
and difficult to eliminate. Moreover, the Monte Carlo pre- 
dictions disfavour a major contribution of long duration 
events. 

Note added in proof. After submission of this work, 
the MEGA collaboration presented their results obtained 
independently from the same data (De Jong et al., 
arXiv:astro-ph/0507286 v2]). Their conclusions are dif- 
ferent from ours. We would like to point out that their 
criticism of our analysis is not relevant because, as stated 
in Section 4.1.2, we choose to only consider for self lensing 
evaluation a population of stars with a standard M/L ra- 
tio, which does not need to account for the total dynamical 
mass nor to reproduce the inner rotation curve. 
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